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Field  Modeling  of  Room  Fires 


William  D.  Davis  Glenn  P.  Forney  John  H.  Klote 


Abstract 

The  application  of  the  Harwell  field  model  FLOW3D  to  model  com- 
partment fires  is  investigated.  Two  experiments  are  chosen  to  model  nu- 
merically. The  first  experiment  is  a single  room  fire  where  the  vertical 
structure  of  the  ceiling  jet  produced  by  the  fire  is  measured  and  the  tem- 
perature response  of  simulated  fusible  links  to  the  ceiling  jet  is  available. 
The  second  experiment  consists  of  three  rooms  with  a fire  in  one  of  the 
rooms.  Temperature  measurements  using  thermocouples  are  available  in 
each  of  the  three  rooms  as  well  as  in  the  corridors  connecting  the  rooms. 
These  two  experiments  provide  an  opportunity  to  investigate  both  two 
and  three  dimensional  field  modeling  of  fires.  It  is  found  that  the  nu- 
merical results  using  the  field  model  are  in  reasonable  agreement  with  the 
experimental  data.  FLOW3D  has  been  enhanced  by  the  addition  of  a sim- 
ple fusible  link  algorithm  previously  used  in  the  zone  fire  model  LAVENT. 
This  algorithm  used  in  conjunction  with  the  field  model  produces  good 
agreement  with  the  measured  fusible  link  temperatures  found  in  the  single 
room  experiment. 


1 Introduction 

The  purpose  of  this  report  is  to  describe  the  use  of  the  Harwell  field  model, 
FLOW3D,  for  modeling  compartment  fires[l] . This  model  was  applied  by  Har- 
well to  simulate  the  fire  at  the  King’s  Cross  underground  station  [2].  In  that 
application  the  field  model,  surprisingly,  predicted  that  the  fire  would  lay  down 
(instead  of  rising  straight  up)  in  the  trench  formed  by  the  walls  and  steps  of  the 
wooden  escalator.  Although  this  prediction  was  met  at  first  with  much  skep- 
ticism, it  was  shown  by  both  the  model  and  by  experiment  that  the  behavior 
of  the  fire  was  due  to  two  related  phenomena.  The  first  effect  also  known  as 
the  ‘Coanda’  effect  was  the  preferential  entrainment  of  air  into  the  plume  from 
below.  The  second  effect  was  the  high  wind  speeds  (14.5  m/s  predicted  by  the 
model  [2])  along  the  base  of  the  escalator.  These  two  phenomena  together  have 
been  referred  to  as  the  “trench  effect” . Similar  insight  may  be  attained  by  ap- 
plying tools  such  as  this  to  other  fire  scenarios  that  have  complex  geometries 
or  flow  fields.  But  first,  further  comparisons  of  the  model’s  predictions  with 
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laboratory  fire  experiments,  which  are  both  controlled  and  well  instrumented, 
are  necessary  to  gain  confidence  in  the  predictive  capability  of  the  model. 

The  experiments  described  in  this  report  serve  as  an  initial  basis  for  verifying 
the  application  of  FL0W3D  for  the  analysis  of  compartment  fires.  The  compar- 
isons between  the  numerical  simulations  and  the  experimental  results  serve  only 
as  a preliminary  verification  of  this  field  model  since  radiation  calculations  were 
not  included  in  the  simulation  and  the  fires  were  modeled  using  experimentally 
measured  heat  release  rates  rather  than  the  available  combustion  models. 

Both  two  and  three  dimensional  approaches  to  modeling  fire  and  smoke  flow 
were  investigated.  Two  experiments  were  chosen  to  compare  with  the  field 
model.  The  first  experiment  [3]  was  selected  because  the  structure  of  the  ceiling 
jet  was  measured  in  the  fire  room.  The  flow  fields  could  be  approximated  using 
a two  dimensional  grid  since  the  experiment  consisted  of  a single  room  with  a 
wide  door.  The  numerical  grid  was  chosen  fine  enough  to  resolve  the  ceiling  jet 
structure.  The  experiment  also  had  a series  of  simulated  fusible  links  located 
at  a variety  of  distances  from  the  ceiling  which  provided  a basis  for  comparing 
a simple  model  of  fusible  link  heating  with  experimental  data. 

The  second  experiment  [4]  was  chosen  because  its  instrumentation  and  ge- 
ometry provided  data  to  illustrate  three  dimensional  variations  in  the  flow  field. 
It  was  a three  room  experiment  that  had  been  used  earlier  to  measure  the  cal- 
culation accuracy  of  the  zone  fire  model  FAST  [5].  Here,  the  purpose  was  to 
investigate  the  effect  of  using  a coarse  grid  in  a field  modeling  calculation.  Un- 
like the  first  experiment,  no  effort  was  made  to  resolve  the  ceiling  jet  and  fusible 
links  were  not  included  in  the  calculation.  The  ability  to  use  a coarse  grid  in 
modeling  a fire  experiment  is  important  due  to  the  limitation  of  computer  speed 
which  places  a practical  limit  on  the  number  of  nodes  used  in  a grid,  particularly 
for  a three  dimensional  calculation. 

1.1  The  Field  Model 

The  fire  experiments  discussed  in  this  report  were  simulated  using  release  2.3.2 
of  FLOW3D  [1].  This  model  has  been  developed  over  a period  of  roughly  a 
decade  at  Harwell  Laboratories.  The  field  model  solves  the  three  dimensional 
form  of  the  equations  for  conservation  of  mass,  momentum,  and  energy  and  the 
k — € model  is  used  to  model  turbulence.  The  required  physical  parameters 
for  these  equations  include  the  fluid  density,  pressure,  specific  heat  at  constant 
pressure,  acceleration  of  gravity,  thermal  conductivity  and  molecular  viscosity. 
The  structure  of  the  field  modeling  computer  package  includes  pre-processing 
and  post-processing  software  as  well  as  the  field  model  itself.  The  interaction 
between  these  three  portions  of  this  package  and  the  files  they  both  produce 
and  require  is  illustrated  in  figure  1. 
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Figure  1:  FL0W3D  Program  and  Data  File  Dependency  Diagram 
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1.1.1  Pre-Processing  Software:  GRIDGEN 

The  pre-processing  software  consists  of  a numerical  grid  generation  package 
named  GRIDGEN  which  is  used  to  generate  body-fitted  coordinates.  This  pack- 
age was  not  used  to  generate  a grid  for  the  experiments  modeled  in  this  report 
since  the  command  file  used  to  run  the  field  model  has  an  option  which  lets 
the  user  create  simple  grids.  This  facility  was  more  than  adequate  to  create  the 
simple  grids  used  for  these  experiments.  Complex  grids  must  be  generated  using 
the  numerical  grid  generator.  The  numerical  grid  generator  uses  a free  format 
input  file  which  allows  the  user  to  define  the  grid  generation  parameters  using 
keywords.  Special  curves  may  be  generated  by  using  a FORTRAN  program  to 
generate  the  necessary  node  pattern  which  is  then  used  by  the  grid  generator 
to  construct  the  required  grid. 

GRIDGEN  creates  a file  named  grid. data  that  is  used  by  the  field  model  to 
determine  the  location  of  the  grid  points. 

1.1.2  Processing  Software:  FLOW3D 

To  use  the  field  model,  a free-form  command  file  must  be  created  which  uses 
keywords  to  set  the  program  options.  Command  files  used  for  modeling  both 
experiments  are  given  in  appendix  A.  In  addition  to  the  command  file,  flex- 
ibility in  setting  program  options  is  provided  by  the  “FORTRAN  frontend”. 
This  frontend  consists  of  a series  of  FORTRAN  subroutines  which  the  user  may 
modify  to  change  physical  parameters,  add  physics,  modify  initializations  and 
boundary  conditions,  or  change  the  format  of  the  program  output.  The  sub- 
routines are  re-compiled  and  re-loaded  to  form  a new  executable1  each  time  the 
field  model  is  run  giving  the  user  the  flexibility  to  customize  its  use.  The  fusible 
link  activation  algorithm  was  added  via  the  FORTRAN  frontend.  The  names 
and  a brief  description  of  each  subroutine  used  in  the  FORTRAN  frontend  for 
both  experiments  are  given  in  appendix  B. 

The  field  model  produces  several  output  files  designated:  dump,  out,  srf  and 
topi.  The  dump  file  contains  the  simulation  results  in  a binary  format.  The 
out  file  is  examined  by  the  user  to  determine  whether  the  run  was  successful,  it 
whether  mass  residuals  were  sufficiently  small.  The  srf  file  contains  information 
used  by  SURF3D.  The  topi  file  contains  topology  information  describing  grid 
connections.  If  the  grid  was  generated  by  the  field  model  then  it  creates  a file 
named  grid. data  containing  grid  point  locations.  These  files  (except  for  the  out 
file)  are  used  by  the  post-processing  software  described  next. 

1.1.3  Post-Processing  Software:  OUTPROC,  TRAN  and  SURF3D 

The  post-processing  software  includes  three  output  graphics  programs:  OUT- 
PROC, TRAN,  and  SURF3D.  The  first  program,  OUTPROC,  produces  contour 

1 An  executable  is  a file  containing  machine  language  instructions  executed  by  the  computer. 
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plots,  velocity  vector  plots,  and  particle  tracks  in  arbitrary  planes  throughout 
the  computational  space.  OUTPROC  uses  a free  format  input  file  which  allows 
the  user  to  define  the  characteristics  of  the  generated  plots  using  keywords. 
An  example  of  the  input  file  is  given  in  appendix  A for  the  two  dimensional 
calculation. 

The  second  output  graphics  program,  TRAN,  provides  a two  dimensional 
animation  of  the  time  evolution  of  a designated  variable  on  the  grid.  Each  node 
point  on  the  grid  is  given  a color  representing  the  magnitude  of  the  variable  at 
that  node.  The  program  displays  the  variable  values  at  each  timestep. 

The  subroutines  necessary  to  generate  a TRAN  file  are  listed  in  appendix 
B.  The  files  produced  by  these  subroutines  are  T.data,  U.data,  P.data  etc  for 
temperature,  velocity  in  the  x direction,  pressure.  A file  may  be  produced 
for  each  of  the  solution  variables.  Up  to  four  variables  may  be  simultaneously 
tracked  using  this  software  with  the  choice  of  variables  specified  within  the 
TRAN  subroutines.  There  is  no  user  interface  for  the  TRAN  files,  hence  any 
modification  to  the  output  must  be  done  by  the  user  directly  to  the  TRAN 
files.  Before  TRAN  is  run,  a program  called  CONVERT  is  used  to  convert  files 
produced  by  the  field  model  (T.data  etc ) into  a file  named  trans.data.  TRAN 
uses  this  file  to  produce  an  animated  display  for  one  of  the  solution  variables. 
At  present,  TRAN  can  only  be  used  for  two  dimensional  calculations. 

The  third  output  graphics  program,  SURF3D,  provides  a three  dimensional 
visualization  of  the  surface  boundaries  (hence  the  name  SURF3D)  at  one  par- 
ticular time  for  the  grid  used  in  the  calculation.  Several  variables  such  as  tem- 
perature or  pressure  may  be  displayed  on  the  grid  boundaries  using  a color  rep- 
resentation for  the  variable’s  magnitude.  The  grid  may  be  expanded,  rotated, 
or  translated  and  the  user  has  the  option  to  move  within  the  grid  boundaries 
to  observe  the  internal  structure  of  the  grid  surfaces.  Before  SURF3D  is  run,  a 
program  named  PRESURF  is  used  to  convert  files  generated  by  the  field  model 
into  a file  used  by  SURF3D  to  generate  the  three  dimensional  display. 

1.1.4  General  Field  Model  Setup  Information 

The  physical  setup  for  each  simulation  was  as  follows.  It  was  assumed  that  the 
fluid  was  air  and  that  it  was  fully  compressible.  Turbulence  was  modeled  using 
the  k — e model  and  included  the  default  treatment  of  the  boundary  layer  at  the 
walls.  The  walls  were  assumed  to  be  adiabatic  except  in  one  three  dimensional 
calculation  where  heat  conduction  in  the  ceiling  was  included.  Radiation  effects 
were  not  included  in  the  calculation  except  to  reduce  the  experimental  heat 
release  rates  by  thirty-five  percent  to  simulate  radiation  loss  to  the  walls  in  the 
fire  room. 

A uniform  grid  was  used  to  model  each  experiment.  All  terms  in  the  equa- 
tions except  for  the  advection  terms  were  discretized  in  space  using  second- 
order  centered  differencing.  For  the  advection  terms,  a number  of  discretization 
schemes  were  available.  For  this  work,  the  upwind  differencing  scheme  was  used. 
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The  equations  were  advanced  in  time  using  a fully  implicit  backward  difference 
procedure. 

1.2  Field  Model  - Zone  Model  Comparisons 

While  this  report  is  primarily  concerned  with  field  modeling,  it  is  important  to 
describe  what  a zone  model  is  in  relation  to  a field  model  since  comparisons 
with  zone  model  calculations  will  be  made.  As  mentioned  earlier,  a field  model 
solves  a basic  set  of  equations  for  conservation  of  mass,  momentum,  and  energy 
on  a user  defined  grid  with  the  coarseness  of  the  grid  determining  the  physical 
scale  lengths  which  can  be  modeled.  The  mass,  enthalpy,  flow  speeds,  turbulent 
energy  etc  are  calculated  for  each  cell  or  control  volume  in  the  grid.  Techniques 
for  modeling  physical  phenomena  whose  natural  scale  sizes  are  smaller  than 
the  grid  size  are  used  in  both  field  and  zone  modeling.  Models  using  these 
techniques  are  sometimes  called  sub-grid  models.  Plume  and  vent  flow  models 
are  examples  of  this  in  zone  models  and  turbulence  models  are  examples  of  this 
in  field  models. 

A zone  model  typically  models  an  enclosure  by  dividing  the  volume  of  the 
enclosure  into  two  layers  or  zones,  designated  the  upper  and  lower  layer.  A 
field  model  divides  the  enclosure  into  many  control  volumes.  As  with  field 
models,  the  temperature  and  density  are  assumed  to  be  uniform  throughout  the 
control  volume  (a  layer  in  the  zone  modeling  case,  an  arbitrary  grid  cell  in  the 
field  modeling  case)  at  a given  time.  The  conservation  equations  for  mass  and 
energy  are  solved  for  each  layer.  The  conservation  of  momentum  equation  is 
usually  neglected  in  a zone  model. 

Perhaps  the  most  notable  difference  between  a field  and  a zone  model  (be- 
sides the  number  of  control  volumes  used)  is  the  method  used  to  model  the 
exchange  of  mass  and  enthalpy  between  control  volumes.  A field  model  solves 
for  these  flows  directly  using  the  conservation  laws  while  a zone  model  uses  var- 
ious sub-grid  models  (plume  and  vent  flow  for  example)  to  model  the  transfer 
of  mass  and  enthalpy  from  the  fire  to  the  layers,  the  entrainment  of  lower  layer 
gas  into  the  upper  layer  and  the  exchange  of  mass  and  enthalpy  between  con- 
nected spaces  through  vents  etc.  Transient  effects,  such  as  the  time  required  to 
establish  the  plume  or  an  upper  layer,  are  usually  ignored  in  zone  models. 

A zone  model  typically  neglects  the  hydrostatic  terms  in  the  equation  of 
state  (the  ideal  gas  law)  which  relates  density,  temperature  and  pressure  since 
the  contribution  of  the  hydrostatic  term  is  small.  However,  these  terms  are  not 
neglected  when  calculating  vent  flows.  These  assumptions  concerning  hydros 
static  terms  are  discussed  in  [6]. 

The  effect  of  gas  flow  in  either  layer  is  ignored  with  the  exception  that  some 
zone  models  contain  a set  of  ceiling  jet  equations  which  represent  the  flow  of  the 
hot  plume  gases  across  the  ceiling.  Transient  effects  of  setting  up  the  ceiling  jet 
are  ignored  in  these  models. 

The  flow  of  the  layer  gases  through  vents  to  other  enclosures  is  calculated 
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using  Bernoulli’s  equation  with  a vent  coefficient  where  the  hydrostatic  pressure 
difference  across  the  vent  as  a function  of  height  is  obtained  from  the  layer 
interface  height  and  the  gas  density  in  each  layer.  The  inflowing  or  outflowing 
mass  and  enthalpy  are  then  placed  in  one  of  the  layers  depending  on  the  relative 
temperature  of  the  flowing  gas  and  the  layers. 

The  pressure  field  in  a field  model  is  calculated  at  each  grid  cell  from  the 
conservation  equations.  The  pressure  field  in  a zone  model  is  derived  from  a 
calculated  reference  pressure  (usually  the  pressure  at  the  base  of  the  lower  layer), 
the  layer  interface  height  and  densities  of  both  layers.  Due  to  neglecting  the 
momentum  equation,  the  pressure  at  a given  height  in  a zone  model  is  constant. 
This  need  not  be  true  in  a field  model.  As  a result,  a field  model  can  model 
regimes  with  complex  flow  patterns  whereas  a zone  model  cannot. 

As  mass  and  enthalpy  are  added  to  or  subtracted  from  the  upper  and  lower 
layers  in  a zone  model,  the  interface  height  of  the  two  layers  is  adjusted  to 
maintain  pressure  balance  between  the  layers. 

Examples  of  zone  models  which  contain  fairly  sophisticated  physics  include 
FAST  [5],  FIRST  [7],  and  LAVENT  [8].  Zone  models  are  most  useful  when 
the  fire  room  and  adjacent  rooms  have  a geometry  that  is  a simple  rectangular 
parallelepiped  and  the  flows  introduced  into  the  rooms  by  the  fire  are  not  com- 
plicated. Field  models  will  allow  the  simulation  of  more  complicated  geometries 
and  can  be  used  in  cases  where  the  flow  fields  are  complicated.  Zone  models  re- 
quire much  less  computer  time  to  simulate  a fire  than  field  models.  Hence  zone 
models  are  useful  for  simulating  fire  scenarios  consisting  of  relatively  simple  flow 
fields  and  geometry.  Use  of  a field  model  may  be  required  if  the  fire  simulation 
requires  either  physics  or  geometries  that  zone  models  can’t  emulate. 

2 One  Room  Experiment 

2.1  Experimental  and  Field  Model  Setup 

The  one-room  experiment  described  in  [3]  was  chosen  to  compare  with  the  field 
model  for  two  reasons.  First,  with  the  fire  located  at  the  back  of  the  room  and 
with  the  door  being  both  wide  and  in  the  front  of  the  room,  the  experiment  could 
be  satisfactorily  modeled  in  two  dimensions.  Second,  the  experiment  measured 
the  ceiling  jet  temperature  and  velocity  profiles  beneath  the  ceiling  and  since 
the  room  could  be  modeled  using  two  dimensions  a fine  grid  could  be  used  to 
resolve  the  ceiling  jet.  The  room  dimensions  are  given  in  table  1. 

This  calculation  could  also  be  used  to  test  the  capability  of  the  outlet  bound- 
ary condition  used  by  the  field  model  to  approximate  the  flow  through  an  outlet. 
The  assumptions  made  and  equations  used  for  the  outlet  boundary  condition 
are  given  in  appendix  C.  Two  runs  were  made,  one  with  the  fire  room  connected 
directly  to  the  outside  through  an  outlet  whose  size  was  equal  to  the  height  of 
the  door  and  a second  where  a second  room  of  size  equal  to  the  first  was  used 


7 


Table  1:  Room  Dimensions  for  the  One- Room  Experiment 


Location 

Dimension  (m) 

Room 

2.44  W x 3.66  L x 2.44  H 

Room  doorway 

.76  W x 2.03  H 

in  order  to  move  the  outlet  boundary  a substantial  distance  away  from  the  fire 
room.  The  outlet  boundary  for  the  two  room  situation  was  on  the  far  side  of  the 
second  room  and  extended  the  entire  vertical  height  of  the  wall.  The  door  to 
the  fire  room  then  opened  into  the  second  room  and  was  modeled  using  internal 
solid  nodes  to  approximate  the  transom  depth  of  the  door.  The  second  room 
was  added  to  avoid  modeling  the  door  flow  using  a mathematical  boundary 
condition. 

Near  the  center  of  the  room,  a “tree”  of  20  mil  diameter  chromel-alumel, 
bare-bead  thermocouples  located  at  the  ceiling  and  at  0.152  m intervals  below 
the  ceiling  were  used  to  obtain  the  temperature  in  the  enclosure  as  a function  of 
height.  Also  near  the  center  of  the  enclosure,  a series  of  simulated  sprinkler  link 
disks,  each  having  an  RTI2  of  25.7  (ms)  5,  diameter  of  0.0095  m,  and  thickness  of 
0.00079  m were  positioned  at  0.01,  0.05,  0.17,  0.22,  and  0.33  m below  the  ceiling. 
It  should  be  noted  that  reference  [3]  reports  an  RTI  of  25  (ms)  5 for  the  sprinkler 
link  disks  but  after  rechecking  the  data,  the  RTI  of  the  disks  ranged  from  25.5 
to  25.8  (ms)  5 [9].  The  total  heat  release  of  the  combustibles  as  a function 
of  time  was  deduced  from  oxygen  consumption  calorimetry  measurements  in  a 
collection  head  outside  the  door. 

Three  different  heat  release  rate  curves  [10]  have  been  used  to  describe  the 
fire  which  was  a waste-paper  basket  initiated  office  furniture  fire  in  one  corner 
of  the  room.  The  linear  heat  release  curve  was  used  based  on  a successful 
simulation  of  the  fire  using  the  zone  computer  model  LAVENT  [10].  The  fire 
volume  was  estimated  using  the  heat  release  rate  for  a fully  involved  fire  of  wood 
pallets  [11].  The  volume  of  the  fire  was  kept  constant  throughout  the  simulation 
but  the  heat  release  rate  per  cell  was  varied  to  mimic  the  total  heat  release  of 
the  fire. 

It  was  assumed  that  35%  of  the  heat  release  was  lost  to  radiation  to  the 
walls.  Since  the  walls  were  considered  adiabatic,  that  meant  that  only  65% 
of  the  heat  release  rate  went  directly  to  heating  the  gas.  This  ratio  was  kept 
constant  throughout  the  simulation  of  the  fire. 

2 Response  time  index 
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Figure  2:  Setup  for  the  One-Room  Experiment 


The  fusible  links  were  simulated  using  the  simplified  model  of  Heskestad 
and  Smith  [12].  Here,  the  heating  of  the  links  requires  the  calculation  of  both 
the  flow  velocity  and  temperature  at  the  link  location.  The  equation  used  to 
compute  the  link  heating  is  given  by 

dT\  _ (T  — T\)\fV 
~dt  ~ RTI 

where  T\  is  the  link  temperature,  T is  the  gas  temperature,  V is  the  gas  velocity, 
and  RTI  is  the  response  time  index  of  the  link.  This  model  neglects  the  effects 
of  radiation  to  the  links  or  conduction  from  the  links  to  the  ceiling. 

Since  the  fusible  links  were  located  in  the  center  of  the  room,  the  diagonal 
distance  from  the  fire  center  to  the  link  position  was  used  to  model  the  flow 
distance.  No  effort  was  made  to  model  the  links  near  the  door  because  of  the 
expected  three  dimensional  distortion  of  the  flow  field  near  at  the  door. 

The  experiment  was  modeled  using  a grid  of  42  x 42  in  the  xy  plane  with  the 
y dimension  oriented  vertically  upward  (see  figure  2).  This  grid  actually  models 
two  rooms  of  equal  size  with  the  fire  located  on  the  far  right  section  of  the 
first  room  and  the  vent  occupying  the  entire  wall  of  the  left  side  of  the  second 
room.  The  additional  room  was  added  to  eliminate  the  observed  influence  that 
the  vent  boundary  conditions  introduced  into  the  calculations  when  a grid  of 
22  x 42  was  used  to  simulate  the  experiment.  These  boundary  effects  will  be 
discussed  later  in  the  report. 

Heat  conduction  to  the  ceiling  was  ignored.  This  latter  approximation  was 
reasonable  since  the  walls  were  made  of  marinite  which  is  a good  insulator. 
The  heat  flow  through  the  walls  during  the  200  seconds  of  the  experiment  was 
therefore  of  marginal  importance.  Radiation  to  the  walls  and  the  gas  was  also 
ignored  as  mentioned  earlier. 
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Figure  3:  Ceiling  Jet  Temperatures  at  100  Seconds 


The  problem  was  run  using  a timestep  of  .5  seconds.  Larger  timesteps  pro- 
duced poor  convergence  for  the  solver  with  resulting  large  residuals.  Details  of 
the  command  file  are  given  in  appendix  A. 

2.2  Comparison  of  Experimental,  Zone  Modeling  and  Field 
Modeling  Results 

The  ceiling  jet  temperatures  predicted  by  FLOW3D  and  LAVENT  are  compared 
with  the  experimental  values  as  shown  in  figure  3 and  4.  As  can  be  observed, 
both  the  field  and  zone  models  give  a good  representation  of  the  temperature  of 
the  jet  at  100  and  150  seconds.  A comparison  of  the  fusible  link  temperatures 
as  calculated  by  these  two  programs  versus  experiment  is  given  in  figures  5 and 
6 for  100  and  150  seconds.  Again,  both  models  provide  a good  representation 
of  the  fusible  link  heating.  It  should  be  noted  that  both  programs  use  the  same 
physical  model  for  fusible  link  heating.  However,  the  ceiling  jet  calculation  in 
LAVENT  is  modeled  using  a zone  model  approximation  [13]. 

The  zone  model  calculates  an  average  layer  height  which  is  uniform  over 
the  entire  room.  The  field  model  on  the  other  hand  calculates  the  flow  field 
dynamics  and  hence  the  predicted  hot  smoke  layer  depth  varies  from  one  place 
in  the  room  to  another.  In  particular,  it  is  pulled  upward  and  entrained  in 
the  vicinity  of  the  plume  as  shown  in  the  temperature  contour  plot  in  figure 
7.  Part  of  the  hot  smoke  layer  is  also  deflected  by  the  top  of  the  door  and 
part  of  the  deflected  layer  is  directed  back  toward  the  plume.  This  flow  of 
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Figure  4:  Ceiling  Jet  Temperatures  at  150  Seconds 


Figure  5:  Fusible  Link  Temperatures  at  100  Seconds 
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T TEMPERATURE  CONTOURS 


Figure  7:  Temperature  Contour  Plot 
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Figure  8:  Velocity  Vector  Plot 


the  hot  smoke  gases  back  toward  the  plume  is  best  observed  using  the  output 
package  “TRAN”  which  provides  a time  dependent  view  of  the  development  of 
the  gas  temperature  during  the  computation.  The  time  independent  plot  of  the 
velocity  flow  shown  in  figure  8 does  not  show  the  flow  back  toward  the  plume 
in  a convincing  manner.  The  arrow  size  in  figure  8 scales  with  the  flow  velocity 
with  the  larger  arrows  representing  greater  speeds.  The  direction  that  the  arrow 
points  gives  the  direction  of  the  vector  flow.  The  mixing  of  the  hot  outflowing 
gas  with  the  cool  inflowing  gas  is  a feature  which  is  only  modeled  in  the  most 
sophisticated  of  the  zone  models. 

The  effect  of  the  boundary  condition  was  studied  by  comparing  the  following 
two  numerical  experiments.  First,  the  door  was  modeled  using  the  boundary 
condition  in  a one  room  simulation  of  the  experiment.  Second,  the  door  was 
modeled  using  a solid  grid  to  represent  the  door  soffit  with  the  boundary  moved 
to  the  far  end  of  a second  room.  The  values  for  the  horizontal  and  vertical 
flow  velocity  along  the  door  at  the  grid  points  just  inside  the  door  were  then 
compared  for  the  two  calculations.  A plot  of  the  horizontal  and  vertical  flow 
velocities  along  the  door  for  the  grid  just  inside  the  door  is  given  for  the  one 
and  two  room  calculations  (figures  9 and  10)  at  200  seconds.  Grid  number 
2 represents  the  control  volume  just  above  the  floor  with  the  top  of  the  door 
beginning  at  grid  number  36.  As  can  be  seen,  there  are  some  small  differences 
in  the  x and  y velocity  profiles  with  the  x flow  velocities  showing  the  largest 
deviations  between  the  two  methods  of  modeling  the  door.  The  largest  contrast 
between  the  two  cases  can  be  shown  if  a plot  of  temperature,  figure  11,  along  the 
grid  just  inside  the  door  is  examined.  The  outlet  boundary  condition  (one  room 
calculation)  ignores  any  heating  of  inflowing  gas  by  mixing  with  the  outflowing 
gas  and  brings  in  the  gas  at  the  boundary  reference  temperature.  The  net  result 
is  that  the  temperature  of  the  inflowing  gas  is  underestimated  using  the  outlet 
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Figure  9:  Horizontal  Vent  Velocities  at  200  Seconds 


boundary  condition  which  then  leads  to  a lower  room  temperature  for  the  fire. 
The  outlet  boundary  condition  should  only  be  used  far  from  the  fire  source 
where  temperatures  are  low  and  heating  from  the  mixing  of  gas  at  the  outlet 
would  be  unimportant  to  the  overall  fire  simulation. 

The  numerical  results  also  allow  a comparison  of  the  zone  model  approxima- 
tion for  the  vent  flows  with  that  calculated  by  the  field  model.  The  zone  model 
approximation  assumes  that  the  flow  through  the  vent  is  a function  of  the  hydro- 
static pressure  difference  across  the  vent  and  a vent  coefficient.  The  hydrostatic 
pressure  in  the  zone  model  approximation  will  be  a function  of  layer  tempera- 
ture and  height  within  the  layer.  Since  the  layer  temperature  is  assumed  to  be 
a constant  in  the  zone  model,  the  pressure  within  a layer  will  vary  linearly  with 
height  within  the  layer.  Figure  12  gives  a plot  of  the  x direction  flow  velocity 
versus  flow  velocity  if  only  a hydrostatic  head  is  used  where  the  pressure  beyond 
the  door  is  assumed  to  remain  at  ambient.  The  layer  height  was  obtained  from 
the  field  modeling  calculations  of  the  two  room  experiment  using  the  method 
proposed  by  Cooper  [14]  and  the  layer  temperatures  were  then  calculated  by 
averaging  the  temperature  in  each  layer.  The  shape  of  the  inflow  velocities  as 
calculated  by  the  zone  model  approximation  compares  favorably  with  the  field 
model  results.  If  the  layer  height  had  been  fixed  at  a lower  level,  the  comparison 
of  the  shapes  would  be  excellent. 

The  deviation  between  the  field  model  and  zone  model  calculations  at  node 
35  and  above  occurs  because  of  presence  of  the  door  soffit  which  was  included  in 
the  field  model  calculation  but  was  ignored  in  the  zone  model  calculation.  The 
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Figure  11:  One  and  Two  Room  Temperature  Profile  in  the  Vent 
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Figure  12:  Horizontal  Velocity  Profiles  in  the  Vent 


soffit  acts  as  a solid  barrier  thus  rapidly  reducing  the  x or  horizontal  velocity 
flow  for  the  field  model  calculation.  The  flow  around  the  soffit  also  introduces 
a vertical  component  to  the  flow  which  is  readily  seen  in  figure  10.  The  zone 
model  approximation  assumes  that  the  flow  velocities  are  only  horizontal  with 
no  vertical  component. 

One  of  the  major  differences  between  the  field  model  calculation  and  the 
zone  model  approximation  is  that  the  zone  model  gives  the  location  of  the 
maximum  inflow  velocity  to  be  at  the  bottom  of  the  door  while  the  field  model 
calculation  shows  that  the  maximum  inflow  velocity  occurs  above  the  bottom 
of  the  door.  The  reason  for  the  difference  is  that  the  zone  model  only  includes 
pressure  differences  to  determine  flow  velocities  while  the  field  model  includes 
not  only  pressure  differences  but  also  momentum  and  viscous  effects. 

A comparison  of  the  predicted  zone  model  temperature  profile  compared 
with  the  calculated  two  room  temperature  profile  is  shown  in  figure  13.  The 
zone  model  temperatures  are  deduced  by  finding  an  average  temperature  in  each 
layer  after  the  layer  height  has  been  calculated.  Calculating  the  layer  height 
properly  is  critical  in  determining  whether  the  enthalpy  flow  out  the  vent  will 
be  successfully  approximated  by  a zone  model. 

2.3  Computational  Accuracy 

A study  of  the  effects  of  computational  accuracy  was  performed  for  this  case.  To 
introduce  the  computational  methods  used  in  the  field  model,  the  conservation 
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Figure  13:  Zone  and  Field  Model  Temperature  Profiles  in  the  Vent 


or  Navier-Stokes  equations  are  discretized  to  form  a set  of  non-linear  finite 
difference  equations.  This  is  done  by  replacing  each  differential  operator  with 
a finite  difference  analog.  Central,  backward,  and  up-wind  difference  formulas 
are  three  discretizations  used  by  the  field  model.  The  problem  then  is  to  find  at 
various  times  the  mass,  momentum  (for  each  direction),  enthalpy  and  if  a k — e 
turbulence  model  is  used  a turbulent  kinetic  energy  and  a turbulent  dissipation 
rate  that  satisfies  the  non-linear  difference  equations  for  each  control  volume  or 
node.  This  set  of  non-linear  equations  is  solved  iteratively  by  making  a series 
of  guesses.  Each  iteration  of  the  non-linear  problem  is  called  an  outer  iteration. 
The  non-linear  difference  equations  are  linearized.  A linear  system  must  then 
be  solved  for  each  solution  variable,  i.e.  velocity,  enthalpy  etc.  Again,  an 
iterative  procedure  is  used.  Each  iteration  for  these  linear  problems  is  called  an 
inner  iteration.  Iterative  methods  are  necessary  for  solving  the  outer  non-linear 
difference  equations  and  the  inner  linear  systems  due  to  the  large  number  of 
equations  that  must  be  solved.  Direct  methods  such  as  Gaussian  elimination 
are  not  feasible  except  for  trivially  small  problems. 

The  number  of  inner  and  outer  iterations  affect  both  the  accuracy  of  the 
results  and  the  amount  of  computer  time  required  to  obtain  them.  The  residual 
mass  at  each  time  step  is  the  sum  of  all  the  net  mass  flow  between  control 
volumes  in  the  system.  If  there  are  no  sources  or  sinks  then  the  sum  of  all 
mass  exchanges  between  control  volumes  must  be  zero.  It  is  essential  that  the 
mass  residual  be  sufficiently  small  to  insure  that  the  flow  fields  are  modeled 
accurately.  Figures  14,  15,  and  16  show  the  relationships  between  the  number 
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of  inner  pressure  iterations,  number  of  outer  iterations,  timestep  size,  computer 
time  or  time  to  perform  the  calculation  and  the  size  of  the  mass  residuals.  There 
is  a tradeoff  between  accuracy  and  efficiency. 

Figure  14  shows  the  importance  of  solving  the  pressure  equation  accurately. 
Increasing  the  inner  pressure  iterations  from  20  to  50  results  in  a three  order 
of  magnitude  reduction  in  mass  residuals.  Figure  15  shows  that  as  the  number 
of  outer  iterations  increases  the  mass  residuals  decrease  and  the  computer  time 
increases.  Figure  16  illustrates  the  effect  of  time  step  size  on  computer  time  and 
mass  residual  size.  As  the  time  step  size  increases  the  residual  size  increases  and 
the  computer  time  decreases.  The  computer  time  required  to  solve  the  problem 
is  inversely  proportional  to  the  time  step  size. 

The  computational  experiments  were  performed  on  a Convex  Cl 20  and  a 
Silicon  Graphics  Personal  Iris  model  4D35.  Both  computers  used  32  bit  or 
single  precision  arithmetic  for  these  calculations.  Though  the  architectures  of 
these  two  computers  are  significantly  different3  the  computer  times  for  cases 
run  for  this  report  were  comparable  (within  10  per  cent).  Typical  2-d  cases 
described  in  this  report  took  about  seven  hours  to  run.  The  3-d  cases  took  an 
order  of  magnitude  longer,  about  100  hours.  Much  time  was  also  required  to 
set  up  the  required  input  files  and  to  analyze  the  output. 

3 The  Cl 20  has  vectorizing  hardware  that  exploits  operations  involving  matrices,  vectors, 
etc.  The  4D35  has  super-scalar  hardware  consisting  of  the  MIPS  3000/3010  RISC  (reduced 
instruction  set  chip)  chip  set  running  at  36  MHZ 
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Figure  15:  Outer  Iteration  Impact  on 
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Figure  16:  Time  Step-Size  Impact  on  Computer  Time  and  Mass  Residuals 
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Table  2:  Room  Dimensions  for  the  Three-Room  Experiment 


Location  Dimension  (m) 


First  room 

First  room  stub  corridor 
First  room  doorway 
Second  room 

Second  room  exit  doorway 
Third  Room 

Third  room  stub  corridor 
Third  room  doorway 


2.34  W x 2.34  L x 2.16  H 

1.02  W x 1.03  L x 2.00  H 

.81  W x 1.6  H 

2.44  W x 12.19  L x 2.44  H 

.76  W x 2.03  H 

2.24  W x 2.22  L x 2.43  H 

.79  W x .94  L x 2.04  H 

.79  W x 2.04  H 


3 Three  Room  Experiment 

3.1  Experimental  and  Field  Model  Setup 

A set  of  experiments  using  a three  room  configuration  [4]  was  chosen  to  test 
the  three  dimensional  capabilities  of  the  field  model.  This  particular  set  of 
experiments  used  a burner  as  a fire  source  and  therefore  represented  a much 
more  controlled  experimental  environment  than  the  single  room  test  case.  Of 
the  nine  sets  of  data  available  from  the  experiment,  set  four  was  chosen  based 
on  the  quality  of  the  data.  The  room  geometry  is  shown  in  figure  17  with  the 
room  dimensions  given  in  table  2. 

Basically,  the  fire  room,  room  one,  is  connected  to  a long  second  room  or 
corridor  through  a short  corridor  with  a door.  There  is  a door  at  the  far  end 
of  the  second  room  and  a third  room  connects  to  the  second  room  by  a short 
corridor  near  the  far  end  of  the  second  room. 

A set  of  eight  thermocouple  trees  was  used  to  measure  the  temperature 
at  various  locations  in  the  three  rooms.  These  trees  are  depicted  by  crosses 
inside  circles  in  figure  17.  The  locations,  spacing  and  numbers  of  thermocouples 
on  each  tree  are  given  in  table  3 of  reference  3.  For  this  report,  five  of  the 
eight  thermocouple  trees  are  used  to  compare  the  experimental  results  with  the 
numerical  calculation. 

The  choice  of  grid  was  dictated  by  the  detail  of  the  structure  to  be  modeled. 
The  corridors  and  rooms  had  different  ceiling  heights  and  the  door  heights  were 
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Figure  17:  Three  Room  Experimental  Layout  and  Thermocouple  Configuration 
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less  than  the  ceiling  heights.  To  keep  a reasonable  number  of  node  points,  it 
was  decided  that  a grid  of  32  x 16  x 14  = 7168  cells  with  the  z direction  being 
vertical  would  be  used.  This  grid  provided  a spatial  resolution  of  .4m  in  the 
x and  y directions  and  .2  m in  the  z direction.  Since  the  smallest  change  in 
vertical  dimensions  was  .2  m,  all  the  changes  in  height  could  be  modeled  using 
a uniform  grid. 

A diffusion  flame  burner,  centered  at  the  middle  of  the  back  wall  of  the  first 
room  was  the  fire  source.  The  top  side  of  the  burner  had  a .34  x .34  m porous 
ceramic  surface  with  a perimeter  of  13  mm  wide  steel  plate.  The  porous  surface 
was  positioned  .5  m above  the  floor.  A 2.8  kW  pilot  flame  was  attached  to  the 
front  of  the  burner.  The  gas  used  was  a mix  of  natural  gas  and  acetylene  in  a 
77  kW  to  23  kW  ratio  respectively.  The  pilot  flame  was  ignited  first  and  burned 
for  300  seconds  and  then  the  main  burner  was  turned  on  at  a rate  of  100  kW  for 
data  set  4.  The  burner  was  allowed  to  burn  at  a constant  rate  for  the  next  900 
seconds.  The  fire  was  modeled  assuming  that  the  heat  released  to  the  gas  was 
65%  of  the  heat  release  rate  with  the  remaining  35%  being  radiated  to  the  walls. 
Since  the  walls  were  assumed  to  be  adiabatic  in  this  case,  the  radiation  flux  was 
in  effect  removed  from  the  calculation.  The  heat  release  position  for  the  burner 
was  positioned  in  a single  grid  which  was  adjacent  to  the  center  line  and  to  the 
right  when  facing  the  burner  from  the  fire  room  door.  The  slightly  off  center 
line  position  was  necessary  because  an  even  number  of  grids  were  chosen  for  the 
fire  room.  The  pilot  flame  was  located  at  the  same  off  center  line  position  and 
was  located  in  the  grid  adjacent  to  the  floor.  For  both  burners,  it  was  assumed 
that  the  power  output  was  a constant  over  the  interval  that  the  burner  was  on. 

The  experiment  was  first  modeled  by  assuming  that  the  walls  were  adiabatic 
and  that  the  door  at  the  end  of  the  second  room  could  be  sensibly  modeled 
by  using  the  outlet  boundary  condition.  Only  the  first  500  seconds4  of  the 
experiment  was  modeled.  For  times  longer  than  500  seconds,  until  the  gas  burner 
was  turned  off  only  small  changes  in  temperature  occurred  as  the  experiment 
was  approaching  steady  state. 

A second  calculation  was  made  in  order  to  include  energy  losses  to  the  ceiling 
by  making  the  ceiling  conducting.  This  was  done  by  adding  two  additional  layers 
of  grids  to  the  vertical  direction  to  represent  the  solid  thereby  extending  the  grid 
size  to  32  x 16  x 16.  These  solid  grids  had  a spatial  width  of  0.1  m , and  hence 
the  grid  used  for  the  calculation  was  no  longer  uniform.  It  should  be  noted 
that  the  grids  used  to  model  the  ceiling  were  much  thicker  than  the  actual  wall 
material  thus  providing  the  ceiling  with  much  greater  thermal  inertia  than  in 
the  original  experiment.  The  effect  of  this  approximation  means  that  more  heat 
is  lost  to  the  ceiling  than  was  originally  experienced  in  the  experiment.  This 
effect  is  somewhat  offset  by  the  fact  that  the  walls  are  still  considered  adiabatic 
in  the  calculation.  A numerical  instability  was  encountered  when  using  very 
small  grids  to  model  the  ceiling.  Further  work  needs  to  be  done  in  order  to 

4 300  seconds  of  pilot  burner  and  200  seconds  of  the  main  burner 
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Figure  18:  Temperature  Profiles  for  Thermocouple  Tree  1 at  500  Seconds 


properly  model  situations  that  require  substantial  changes  in  grid  size. 

3.2  Comparison  of  Experimental  and  Numerical  Results 

Figures  18  through  22  give  the  temperature  profile  of  the  five  thermocouple  trees 
versus  the  calculated  temperature  profiles  using  an  adiabatic  and  conducting 
ceiling  at  500  seconds.  Both  ceiling  types  seem  to  give  reasonable  temperature 
comparisons  with  the  thermocouple  trees  in  the  fire  room,  the  connecting  cor- 
ridor and  the  first  thermocouple  tree  in  the  second  room  nearest  the  fire  room. 
The  gas  temperatures  in  the  fire  room  both  near  the  floor  and  near  the  ceil- 
ing are  underestimated  in  the  calculation.  This  is  due  to  the  lack  of  radiative 
heating  of  the  ceiling,  walls,  and  floor.  In  particular,  if  the  floor  were  heated 
by  radiation,  the  inflowing  gas  would  be  warmed  by  the  floor  yielding  higher 
gas  temperatures  near  the  floor.  There  should  also  be  a similar  increase  in  gas 
temperature  near  the  ceiling  due  to  radiation  heating  of  the  walls  and  ceiling. 

For  the  thermocouple  trees  at  the  far  end  of  room  two  and  in  the  corner  of 
room  three  (tree  5 and  tree  8),  the  calculation  overestimates  the  temperature 
of  the  gas.  Part  of  the  reason  for  overestimating  the  temperature  would  come 
from  the  absence  of  heat  loss  to  the  walls.  The  effect  of  allowing  the  ceiling  to 
remove  heat  is  clearly  demonstrated  by  noting  the  decrease  in  gas  temperature 
that  occurred  when  the  ceiling  was  made  conducting  rather  than  adiabatic  as 
shown  in  figures  21  and  22.  It  would  be  expected  that  additional  improvement 
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Figure  19:  Temperature  Profiles  for  Thermocouple  Tree  2 at  500  Seconds 


Figure  20:  Temperature  Profiles  for  Thermocouple  Tree  3 at  500  Seconds 
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Figure  21:  Temperature  Profiles  for  Thermocouple  Tree  5 at  500  Seconds 


Figure  22:  Temperature  Profiles  for  Thermocouple  Tree  8 at  500  Seconds 
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in  the  agreement  between  the  calculated  and  experimental  temperatures  would 
be  realized  by  allowing  the  walls  to  be  conducting.  This  modification  was  not 
tried  due  to  the  time  constraints  imposed  on  the  project. 

It  should  be  noted  that  the  exit  door  out  of  room  two  was  modeled  using 
the  outlet  boundary  condition.  From  earlier  discussion,  one  would  expect  this 
boundary  condition  to  lower  the  local  temperature  since  any  heating  of  the 
inflowing  air  by  the  outflowing  air  is  neglected.  In  this  instance,  however,  the 
air  temperatures  are  fairly  low  at  the  door  and  the  door  is  relatively  isolated 
from  the  rest  of  the  experiment.  It  was  felt  therefore  that  during  the  first  200 
seconds  of  the  main  burning  phase  of  the  experiment,  modeling  the  door  by 
the  outlet  boundary  condition  would  only  have  a small  impact  on  the  overall 
calculations  and  the  savings  in  computer  time  by  minimizing  the  number  of 
nodes  would  be  substantial. 

The  impact  on  the  calculation  of  varying  the  number  nodes  was  only  done  in 
the  vertical  direction.  Two  calculations  were  performed.  The  first  calculation 
used  eight  nodes  to  represent  the  vertical  height.  A second  calculation  which  is 
reported  in  detail  in  this  report  used  14  nodes  to  represent  the  vertical  height. 
While  the  additional  detail  provided  by  the  14  node  calculation  was  valuable  in 
comparing  the  experimental  measurements  with  the  calculation,  the  overall  flow 
patterns  and  the  values  of  the  temperatures  were  approximately  reproduced  on 
the  coarser  grid. 

The  calculation  must  still  be  repeated  with  a doubling  of  nodes  in  the  X 
direction  in  order  to  verify  that  the  flow  in  the  corridors  out  of  the  fire  room 
and  into  room  three  are  begin  simulated  properly.  Since  each  doubling  of  nodes 
in  one  of  the  dimensions  increases  the  computation  time  by  a factor  of  four, 
there  is  a practical  limit  to  the  number  of  nodes  that  can  be  used  in  a three 
dimensional  calculation. 


4 Summary 

The  simulation  of  the  single  room  fire  using  FLOW3D  was  successful  in  that  the 
measured  temperatures  for  both  the  thermocouples  and  the  simulated  fusible 
links  agreed  closely  with  the  calculated  values.  The  absence  of  radiation  in  the 
calculations  did  not  seem  to  effect  the  calculation  of  the  ceiling  jet  interaction 
with  the  simulated  fusible  links. 

The  simulation  of  the  three  room  experiment  was  not  quite  as  successful 
in  that  even  though  the  fire  room  was  fairly  well  modeled,  the  temperature 
predictions  for  the  third  room  and  also  for  the  far  end  of  the  second  room  were 
not  as  good.  One  reason  for  the  poorer  simulation  of  these  areas  far  from  the 
fire  could  be  the  use  of  too  few  nodes  in  modeling  the  flow  through  the  doors 
and  down  the  short  corridors.  A second  reason  could  be  the  absence  of  heat 
convection  to  the  walls  since  both  thermocouple  trees  where  the  temperature 
was  significantly  overestimated  were  quite  far  from  the  fire  and  the  dominant 
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loss  of  energy  in  these  areas  would  be  by  convection  to  the  ceiling  and  walls. 

Additional  effort  should  be  placed  in  simulating  the  three  room  experiment 
The  simulation  needs  to  be  repeated  with  more  nodes  defining  the  door  widths 
in  the  connecting  corridors  and  also  in  the  region  of  the  plume.  The  effect  of 
radiation  needs  to  be  included  and  the  walls  should  be  treated  as  conducting 
surfaces  such  that  convective  losses  to  the  walls  would  be  included  in  the  calcu- 
lation. If  a significant  difference  between  the  results  of  the  simulation  and  the 
experimental  measurements  persist,  it  would  then  be  necessary  to  reexamine 
the  experiment  for  some  feature  that  is  not  being  properly  modeled  and  also 
examine  the  k — € turbulence  model  and  the  convective  heat  transfer  calculations 
to  the  walls  and  ceiling. 
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A Field  Model  Command  Files 


A listing  of  the  command  or  input  files  necessary  to  run  field  model,  FLOW3D 
and  the  graphics  package,  OUTPROC,  are  given  below.  Descriptions  of  each 
command  in  these  files  may  be  found  in  the  user’s  manuals[l]. 

A.l  One-Room  Command  Files 

A. 1.1  Field  Model  Command  File  for  the  One-Room  Experiment 


»FL0W3D 

»SET  LIMITS 
#CALC 
NI=42; 

NJ=42 ; 

NK=3 ; 

#ENDGALC 

NUMBER  OF  I POINTS  #NI 
NUMBER  OF  J POINTS  #NJ 
NUMBER  OF  K POINTS  #NK 
END 

» OPTIONS 

CARTESIAN  COORDINATES 
RECTANGULAR  GRID 
TRANSIENT  FLOW 
TEMPERATURE  CALCULATION 
COMPRESSIBLE  FLOW 
TURBULENT  FLOW 
END 

»M0DEL  DATA 
»GAS  CONSTANTS 
TREF  297. 

»PHYSICAL  PROPERTIES 
»STANDARD  FLUID 
AIR 

/*  REFERENCE  TEMPERATURE  20. 

BUOYANCY  REFERENCE  DENSITY  1.23  */ 
REFERENCE  TEMPERATURE  24. 

BUOYANCY  REFERENCE  DENSITY  1.19 
GRAVITY  VECTOR  0.  -9.8  0. 

END 

»PHYSICAL  CONSTANTS 
DEFAULT  SOLID  CONDUCTIVITY  .0259 
END 

»SET  OPTIONS 
BUOYANCY 

TIME  STEPS  400*. 5 
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NSTEP  400 

/*  INITIAL  TIME  100.00  */ 

END 

/*  »SET  INITIAL  GUESS 
RESTORE  VARIABLES 
LAST  DUMP  FILE 
END  */ 

»DIFFERENCING  SCHEME 
ALL  EQUATIONS  ’UPWIND’ 

END 

»SOLVER  DATA 

»PROGRAM  CONTROL 

MAXIMUM  NUMBER  OF  ITERATIONS  30 
OUTPUT  MONITOR  POINT  38  38  2 
MASS  SOURCE  TOLERANCE  l.E-6 
END 

» NUMBER  OF  SWEEPS 
»MAXIMA 
U 10 

V 10 
K 20 

EPSILON  20 
P 50 
END 

»MINIMA 

U 5 

V 5 
P 35 
K 10 

EPSILON  10 
END 

» UNDER  RELAXATION  FACTORS 
U .6 
V .6 
P .6 
TE  .6 
ED  .6 
T .6 
H .6 

/*  FULL  RESIDUAL  PRINTING  */ 

END 

»REDUCTION  FACTORS 
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U .1 
V .1 
P .01 
TE  .1 
ED  .01 
T .1 
H .1 
END 

»CREATE  GRID 

/*  » REGION 

READ  GRID  FROM  DUMPFILE  */ 

»SIMPLE  GRID 
X START  0. 

Y START  0. 

Z START  0. 

DX  40  * .183 
DY  40  * .0625 
DZ  2 * 1.22 
END 

»MODEL  TOPOLOGY 

»PRESSURE  OUTFLOW  BOUNDARIES 
/*  PRESSURE  LOCATIONS (1,9  TO  #NJ,4)  =0 
PRESSURE  (1,9  TO  #NJ,4)  =0.  */ 

PRESSURE  LOCATIONS (1,2  TO  39,2)  =1 
PRESSURE  (1,2  TO  39,2)  =0. 

END 

»SOLID  REGION 

»NON  CONDUCTING  REGION 

SOLID  CONTROL  VOLUMES  (2  TO  41,40  TO  41,2)  = 0 
/*  SOLID  DENSITY  8.81E3 
SPECIFIC  HEAT  1.25E3.  */ 

END 

»NON  CONDUCTING  REGION 

SOLID  CONTROL  VOLUMES  (21,36  TO  39,2)  = 1 
END 


»STOP 


A.  1.2  OUTPROC  Command  File  for  the  One-Room  Experiment 

/**************************************** *********************/ 

/*  EXAMPLE  1:  (EX1.0CMD)  */ 

/*************************************************************/ 
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» OUTPROC 
» SET  LIMITS 

/* - — — */ 

/*  FOR  2D  PLOTS  SET  : */ 

/*  SPACE  DIMENSIONS  - 2 */ 

I*  NO.  OF  ELEMENTS  = MAX(NI*NJ ,NJ*NK,NK*NI)  */ 

/*  NO.  OF  NODES  - MAX(NIP1*NJP1 ,NJP1*NKP1 ,NKP1*NIP1  */ 

/*  NO.  OF  FREEDOMS  - VARIABLES  * NODES  */ 

/* — — — — — —————— — — */ 


SPACE  DIMENSIONS  2 

VARIABLES  7 

ELEMENTS  9000 

NODES  9000 

FREEDOMS  60000 

END 

» SET  OPTIONS 
TWO  DIMENSIONS 
END 

» SET  VARIABLES 

VARIABLE  NAMES  >U’  ’V’  ’W’  >P>  ’TE’  ’ED’  »T» 

END 

» MODEL  DATA 
» CREATE  GRID 

» GENERATE  A USER  GRID 

MASS  CONTROL  VOLUME  VERTICES 

SOLID  SIDE  INTERPOLATION  /*  SET  VARIABLES  TO  0 AT  WALLS  */ 

READ  FROM  STREAM  21  /*  READ  DUMPFILE  FROM  STREEM  21  */ 

SLICE  COORDINATE  3 /*  IE  K PLANES  */ 

PLANE  NUMBER  2 /*  IE  K-  6 */ 

END 

» LOAD  USER  VALUES 
» OUTPUT  DATA 
/*  » SELECT  TIME 
TIME  9. 

END  */ 

» SET  PLOT  OPTIONS 
HEADING  » EXAMPLE  1’ 

USE  LOCAL  EXTREMA 

/*  EXAMPLE  OF  ZOOMING  IN  ON  A REGION 
X LIMITS  OF  ZOOM  1.0  2.0 
Y LIMITS  OF  ZOOM  0.5  0.75 
COMMENTED  OUT  HERE  */ 
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» PLOT  GRID 

CAPTION  ’ GRID’ 

HEADING  ’ EXAMPLE  1’ 

END 

» SET  OUTPUT  POINTS 
» SELECT  ALL  NODES 

» PLOT  VECTORS 
HEADING  ’ ROOM’ 

CAPTION  ’VELOCITY  VECTORS’ 

VARIABLES  ’U’  ’V’ 

PLOT  REFERENCE  VECTOR 
AUTOMATIC  SCALING 
COLOUR  USING  VARIABLE  ’T’ 

NUMBER  OF  INTERVALS  7 

INTERVALS  290.  300.  310.  320.  330.  340.  350. 

COLOURS  ’BLUE’  ’CYAN’  ’GREEN’  ’YELLOW’  ’ORANGE’  ’PINK’  ’RED’ 
END 

» PLOT  CONTOURS 
HEADING  ’ ROOM  ’ 

CAPTION  ’U  VELOCITY  CONTOURS’ 

VARIABLE  ’U’ 

NUMBER  OF  CONTOURS  12 

COLOURS  ’CYAN’  ’GREEN’  ’YELLOW’  ’RED’ 

END 

» PLOT  CONTOURS 
HEADING  ’ROOM’ 

CAPTION  ’V  VELOCITY  CONTOURS’ 

VARIABLE  ’V’ 

NUMBER  OF  CONTOURS  12 

COLOURS  ’CYAN’  ’GREEN’  ’YELLOW’  ’RED’ 

END 

» PLOT  CONTOURS 
HEADING  ’ROOM’ 

CAPTION  ’P  PRESSURE  CONTOURS’ 

VARIABLE  ’P’ 

NUMBER  OF  CONTOURS  15 

COLOURS  ’CYAN’  ’GREEN’  ’YELLOW’  ’ORANGE’  ’RED’ 

END 

» PLOT  CONTOURS 
HEADING  ’ROOM’ 
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CAPTION  ’TE  ’ 

VARIABLE  ’TE’ 

NUMBER  OF  CONTOURS  11 
COLOURS  ’RED’  ’GREEN’  ’YELLOW’ 

END 

» PLOT  CONTOURS 
HEADING  ’ROOM’ 

CAPTION  ’ ED  ’ 

NUMBER  OF  CONTOURS  11 
COLOURS  ’RED’  ’GREEN’  ’YELLOW’ 

END 

» PLOT  CONTOURS 
HEADING  ’ ROOM’ 

CAPTION  ’T  TEMPERATURE  CONTOURS’ 

NUMBER  OF  CONTOURS  7 

CONTOUR  VALUES  300.  310.  320.  325.  330.  335.  340. 

VARIABLE  ’T’ 

COLOURS  ’BLUE’  ’CYAN’  ’GREEN’  ’YELLOW’  ’ORANGE’  ’PINK’  ’RED’ 
END 


» STOP 

» END  OF  DATA 

A. 2 Three-Room  Field  Model  Command  File 

»FL0W3D 

»SET  LIMITS 
#CALC 
NI=32 ; 

NJ=16 ; 

NK=16 ; 

#ENDCALC 

NUMBER  OF  I POINTS  #NI 
NUMBER  OF  J POINTS  #NJ 
NUMBER  OF  K POINTS  #NK 
NUMBER  OF  DIMENSIONS  3 
TOTAL  INTEGER  WORK  SPACE  50000 
TOTAL  REAL  WORK  SPACE  400000 
TOTAL  INTEGER  LIST  SPACE  10000 
TOTAL  REAL  LIST  SPACE  100000 
END 

» OPTIONS 

CARTESIAN  COORDINATES 
RECTANGULAR  GRID 
TRANSIENT  FLOW 
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TEMPERATURE  CALCULATION 
COMPRESSIBLE  FLOW 
TURBULENT  FLOW 
END 

»MODEL  DATA 
»GAS  CONSTANTS 
TREF  293. 

END 

»PHYSICAL  PROPERTIES 
»STANDARD  FLUID 
AIR 

REFERENCE  TEMPERATURE  20. 

BUOYANCY  REFERENCE  DENSITY  1.23 
GRAVITY  VECTOR  0.  0.  -9.8 
END 

»PHYSICAL  CONSTANTS 
DEFAULT  SOLID  CONDUCTIVITY  .12 
END 

»SET  OPTIONS 
BUOYANCY 

TIME  STEPS  400  * .25 
NSTEP  400 

INITIAL  TIME  400. 

END 

»SET  INITIAL  GUESS 
RESTORE  VARIABLES 
LAST  DUMP  FILE 
END 

»DIFFERENCING  SCHEME 
ALL  EQUATIONS  ’UPWIND’ 

END 

»SOLVER  DATA 

»PROGRAM  CONTROL 

MAXIMUM  NUMBER  OF  ITERATIONS  35 
OUTPUT  MONITOR  POINT  30  6 10 
PRESSURE  REFERENCE  POINT  25  13  4 
MASS  SOURCE  TOLERANCE  l.E-6 
END 

» NUMBER  OF  SWEEPS 
»MAXIMA 
U 10 

V 10  * 

K 10 

EPSILON  10 
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P 50 
END 

»MINIMA 
U 2 
V 2 
P 10 
K 2 

EPSILON  2 
END 

» UNDER  RELAXATION  FACTORS 
U .6 
V .6 
P .6 
TE  .6 
ED  .6 
T .6 
H .6 

/*  FULL  RESIDUAL  PRINTING  */ 

END 

» REDUCTION  FACTORS 
U .1 
V .1 
P .01 
TE  .1 
ED  .01 
T .1 
H .1 
END 

»CREATE  GRID 
»REGION 

READ  GRID  FROM  DUMP FILE 
END 

/*  »SIMPLE  GRID 
X START  0. 

Y START  0. 

Z START  0. 

DX  30  * .4 
DY  14  * .4 

DZ  12  * .2  2 * .10 

END  */ 

»MODEL  TOPOLOGY 

» PRESSURE  OUTFLOW  BOUNDARIES 

PRESSURE  LOCATIONS (1,12  TO  13,2  TO  11)  =0 
PRESSURE  (1,12  TO  13,2  TO  11)  =0. 
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END 


»SOLID  REGION 


»NON  CONDUCTING  REGION 

SOLID  CONTROL  VOLUMES  (2  TO  6,2  TO  9,2  TO  15)  = 0 
SOLID  CONTROL  VOLUMES  (13  TO  25,2  TO  7,2  TO  15)  = 1 
SOLID  CONTROL  VOLUMES  (9  TO  29,8  TO  9,2  TO  15)  = 2 
/*  SOLID  CONTROL  VOLUMES  (7  TO  8,8  TO  9,12  TO  15)  = 3 
SOLID  CONTROL  VOLUMES  (30  TO  31,9,12  TO  15)  =4 
SOLID  CONTROL  VOLUMES  (30  TO  31,8,10  TO  15)=5 
SOLID  CONTROL  VOLUMES  (26  TO  31,2  TO  7,13  TO  15)=6  */ 
END 


»CONDUCTING  REGION 
SOLID  CONTROL  VOLUMES 
SOLID  CONTROL  VOLUMES 
SOLID  CONTROL  VOLUMES 
SOLID  CONTROL  VOLUMES 
SOLID  CONTROL  VOLUMES 
SOLID  CONTROL  VOLUMES 
SOLID  DENSITY  740. 
SPECIFIC  HEAT  1.43E3 
END 


(2  TO  31,  10  TO  15,  14  TO  15)  = 0 
(7  TO  8,  8 TO  9,12  TO  15) =1 
(30  TO  31,9,12  TO  15)=2 
(30  TO  31,  8,10  TO  15)=3 
(26  TO  31,2  TO  7,13  TO  15)=4 
(7  TO  12,  2 TO  7,  14  TO  15)=5 


»STOP 

B Field  Model  Subroutines 

The  following  set  of  subroutines  are  used  in  the  simulation  of  the  experiments 
described  in  this  report.  The  subroutines  are  grouped  by  application  and  will 
not  be  described  in  detail  as  adequate  documentation  exists  in  the  code  for  each 
subroutine. 


B.l  TRAN  Subroutines 

Two  subroutines  are  used  to  build  the  data  files  used  by  the  graphics  program 
TRAN.  These  subroutines  are  USRPLT  and  EXTPLT.  The  number  of  variables 
that  may  be  output  at  present  are  four  and  are  chosen  by  specific  inputs  in 
USRPLT. 

B.2  Outlet  Subroutines 

A series  of  six  subroutines  are  used  to  modify  the  conditions  at  the  outlets. 
From  the  standpoint  of  fire  modeling,  if  the  default  outlet  calculation  is  used, 
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the  temperature  of  any  inflowing  gas  through  the  outlet  is  taken  to  be  equal 
to  the  temperature  of  the  last  adjacent  grid  before  the  outlet.  With  this  set  of 
subroutines,  inflow  temperature  is  changed  to  the  reference  temperature  pro- 
vided that  the  calculation  is  done  using  the  compressible  gas  option.  If  the  gas 
is  treated  as  an  incompressible  gas,  the  inflow  reference  temperature  will  be 
0.  These  routines  also  allow  other  changes  to  be  made  to  the  boundary  con- 
ditions but  the  ease  of  making  the  changes  has  not  been  explored.  The  set  of 
subroutines  are  USRPRS,  BCSTMP,  CENTH,  CSCALR,  CT,  and  CTFN. 

B.3  Heat  Source  Subroutine 

Heat  sources,  both  time  dependent  and  time  independent  may  be  inserted  into 
the  subroutine  USRSRC.  This  subroutine  is  called  at  every  calculation  cycle 
and  provides  access  to  all  the  program  variables  and  the  necessary  parameters 
to  control  how  the  program  variables  are  to  be  used.  When  inserting  heat 
sources  into  the  code,  the  most  important  problem  to  overcome  is  to  make  sure 
that  the  volume  heat  release  rate  is  correct  for  the  fire  being  modeled.  Too 
large  or  too  small  of  volume  heat  release  rate  will  yield  incorrect  temperatures 
which  will  then  have  an  effect  on  the  buoyancy  of  the  gas  and  the  accuracy  of 
the  simulation. 

B.4  Fusible  Link  Subroutine 

The  fusible  link  algorithm  was  entered  in  the  subroutine  USRTRN  which  is 
called  at  the  end  of  each  computed  timestep.  The  entire  set  of  program  variables 
are  available  in  this  subroutine.  Any  computation  which  does  not  couple  back 
into  the  main  set  of  program  variables  can  be  performed  as  a side  calculation 
in  this  subroutine. 

B.5  Machine  Dependent  Parameter  Subroutine 

Subroutine  USRMON  was  used  by  Harwell  to  correct  an  underflow  error  in  the 
turbulent  energy  calculation  that  was  uncovered  when  FLOW3D  was  used  on 
the  thirty-two  bit  workstation.  Basically,  certain  variables  are  prevented  from 
taking  on  values  that  are  smaller  than  can  be  handled  by  these  machines.  This 
subroutine  must  be  present  when  FLOW3D  is  run  using  turbulence. 

C Boundary  Conditions 

The  boundary  conditions  at  the  outlet  are  treated  by  assuming  that  the  flow  is 
fully  developed  such  that 


<9$ 

dn 


= n • = 0 
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where  $ represents  one  of  the  variables  such  as  velocity,  turbulent  energy,  tem- 
perature, etc.  For  Cartesian  geometry,  the  boundary  condition  is  simply  satis- 
fied by  copying  the  value  of  phi  from  the  interior  control  volume  adjacent  to  the 
outlet  into  the  neighboring  exterior  control  volumes.  By  using  the  outlet  sub- 
routines described  in  appendix  B,  the  inflowing  enthalpy  or  temperature  may 
be  set  by  the  user  to  values  other  than  the  ones  found  in  the  neighboring  control 
volumes. 
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